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This  Memorandum  daacribaa  a  three-dimensional  numerical  stress-analysis 
computer  program  which  originatad  at  tha  Department  of  Machanical  Enginaaring, 
Imperial  Collage,  University  of  London.  It  is  based  on  the  Boundary  Integral 
Equation  (BIE)  method  and  it  is  used  to  carry  out  stress  analysis  of  three- 
dimensional  component*  of  couples  geometry.  Stress  intensity  factors  can  be 
derived  for  the  case  where  Che  component  is  cracked.  The  body  is  assumed  to 
have  linear  elastic  isotropic  properties.  The  condition  of  plane  strain  is 
assumed  at  all  points  along  the  crack  front  ascept  where  it  intersects  the  free 
surface,  where  plane  stress  assumptions  are  made.  The  program  is  written  in 
standard  FORTRAN  and  it  is  now  operational  on  the  Materials  and  Structures 
departmental  TAX  computer. 


8tep-by-step  instructions  have  been  presented  together  with  examples 
illustrating  the  input  and  output  formats. 


Validation  and  accuracy  checks  of  the  computer  program  have  been  made  by 
solving  two  bench  mark  problems  where  the  BIE  solutions  are  cohered  with  known 
accurate  solutions. 
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i  nmueucnow 

This  Maworandue  describe*  a  coaputar  prograa  BIE3D,  which  la  used  for  solving  thraa- 
diaenaional  alaatoitadca  problaa*  using  tha  Boundary  Intagral  Equation  (BIE)  aathod. 

Ihi*  pro gran  any  ba  oaad  aa  a  nuaerieal  tool  for  tha  atraaa  analysis  of  thraa-dlaanalonal 
craekad  coaponanta  la  llnaar  alaatic  laotroplc  bodlaa. 

The  fundamental*  of  tha  Boundary  Intagral  Equation  aathod  hart  baan  daacrlbad 
aata  naively  alsewbare  (aaa  for  laatanca  Baaarjaa  aad  Buttarfiald1)  and  *111  not  ba 
rapaatad  bar a.  Xha  boundary  surface  of  tha  body  la  divided  Into  aaall  ragioaa  eallad 
alaaantai  la  thla  prograa  BXX3D,  laoparaaatrle  quadrilateral  alaaaata  with  alght-nodaa 
ara  oaad.  Tha  dlatxlbutloa  of  aurfaea  goon* try,  aad  tha  uakaoun  dlaplacaaanta  aad  trac¬ 
tion*,  ara  aaproaaad  la  taraa  of  quadratic  ahapa  functloaa  of  local  coordlaotoa.  After 
the  dlaplacaaanta  aad  atroaaaa  hove  baan  calculatad,  tha  baale  coacapta  of  linear  alaatic 
fracture  aachaalc*  ara  aaad  to  dotaralaa  tha  atraaa  lateaalty  factor*. 

Tba  origlaal  coaputar  prograa  which  *a a  developed  la  tha  Hacbaalcal  Engineering 
Doper taaat  at  Iqarlal  Collage,  Osivaraity  of  London*,  baa  boon  aodiflad  aad  adopted  into 
tha  Matarlala  aad  (true taraa  dopartaaatal  TAX  coaputar.  Thla  Hawn  ran In*  deacrlbaa  ita 
noa  and  lay*  doua  tha  procedure*  for  data  praparatloa  aad  Input. 

la  ardor  to  rail***  tha  nark  load  for  large  problaa*  and  eliadaate  Hunan  error* 
dariag  the  procoaa  of  data  preparation,  a  thraa-dlaanalonal  neah-gaaaratleu  prograa  baa 
boon  developed  la  collaboration  with  Beuthaapton  University3,  aad  adopted  aa  a  data 
preparation  tool!  it  la  projected  bora  aa  aa  option  to  tba  oaera.  Dor  tha  purpoae  of 
obtaining  a  thmo-diaana  local  rapraaantatlen  of  tho  generated  no  ah  a  thraa-dlaanalonal 
graphica  plotting  prograa  haa  boo*  devolopod  at  DAI,  aa  that  the  diacrotlaatlon  data  can 
ba  ohodtad  vl anally  la  taraa  of  olcaaat  topology  oad  nodal  coerdlaataa  before  naa.  Tha 
plotting  prograa  la  writtoa  la  1AIXC  aad  it  tuna  lntaractivaly  on  a  Dewlatt  Packard 
■RM7T  graphic*  tareiaal. 

Tha  orgenlaatloe  of  tba  BXZ3D  prograa  aa  veil  a*  tho  function  of  anch  aubprograc 
la  BIX3D  la  briefly  daacrlbad  la  aoctloo  2. 

la  aoctloo  3,  tho  toekalqu*  of  aodolling  a  problaa,  tha  elaaant  daaign  aad  aub- 
diviaiona  are  daacrlbad  la  detail.  Alao  daacrlbad  la  aactlon  3  ara  tha  guide-llaaa  for 
tha  preparation  of  aa  input  data  fllo.  In  order  to  aako  thla  preparation  aa  alnple  aa 
poaeiblo,  the  coaputar  pragrea*  for  thn  neab-generatlon  oad  graphica  plotting  hove  been 
briefly  daacrlbad  in  aoctloo  3.4  and  lllnatratod  with  srssplo*  la  aoctloo  S.  Tha  full 
description*  of  those  utility  pro grans  la  doacribod  alaoebere3. 

Tha  tachnique  of  datarainiag  the  atraaa  Intensity  factor*  is  contained  in 
•action  4,  using  a  centre •  cracked  apaeinaa  aa  aa  enable.  The  >H  reanlta  are  coopered 
with  the  veleea  obtained  by  Dewnan4  eaiag  tho  boundary  collocation  nsthod.  In  order  to 
estiaat*  the  accuracy  of  the  IXX3D  cenpnter  prograa  the  acme*  oonosutiwtiou  factor*  tor 
a  flat  plat*  with  a  central  aoraal  bole  have  boon  sales latod  for  canpnrison  with  known 
accnree*  eolation*. 
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2  a  3HS£BBBg!fc  cowoth  nocauM  ro»  am  elastic  boot 

BXS3D  is  i  commit  prograa  for  tha  wlyili  of  elastic  behaviour  of  a  three- 
«*Mnaisaal  Mias  «*•  Mthod  of  Boundary  Xacagcal  Equations.  It  la  written  la  stand- 
Md  IMttUI.  The  body  ia  aasuasd  to  have  Isotropic  propart  las.  la  etdar  to  osa  this 
eoapatar  prograa.  tha  surface  of  tha  throa-dlaaoaiooal  body  oust  ba  dlscratlsad  and  divided 
iato  a  aaall  maker  of  eight-node  quadrilateral  alaaaats.  Boa  tha  gaaaatty  aad  tha  vari- 
atioes  of  tha  booadary  displacements  aad  traetiooa  over  thaaa  alaaaats  ara  approaimtaly 
raprasantad  by  quadratic  shapa  faaetioM. 

With  thaaa  antadatim,  tha  gaaaral  latagral  aqoatioa  can  ba  written  aa  a  sat  of 
simultaneous,  llaasr,  algabraic  equations  ahlch  ara  aolrad  sobjact  to  tha  boaodary  coo- 
ditioas  (either  displacaaant  or  strata)  of  tha  preblaa.  Pisplacanaata  aad  traetieas  of 
tha  asrfaea  nodal  points  ara  produced,  froa  ahlch  tha  stroaaaa  at  tha  aarfaeo  aodaa  aad 
tha  displacements  aad  tha  atroaaaa  at  any  required  iatarnal  points  can  ba  cal cu la tod. 
Datsilad  descriptions  of  tha  BX13B  aathod  aad  tha  — — »  fora  lotion  aaad  ia  tha 
analysis  nay  ba  found  ia  Kaf  2. 


Iha  counter  prograa  11X30  has  baaa  written  la  standard  rxOAK  language  and  it  is 
dlrldad  iato  a  as  la  prograa  aad  a  total  of  24  subprogram.  Iha  ordar  of  tha  subprograas 
ara  liatad  balsa  together  ala  a  brief  description  of  each  faactioai 

(1)  man-  Poraa  the  aatris  aad  tha  second  aeOar  by  BIZ  aathod  (is  right-hand  side 

of  BXI  equations) . 

(2)  ■>-  leads  a  asah  of  sight-node  quadrilateral  alaaaats. 

(3)  W8HOOT-  Writes  oat  tha  Mb  data. 

(4)  ABBAT-  Pats  alasMat  node  ambers  ead  nodal  coordinates  ia  arrays. 

(3)  TBAM-  Treasforae  cylindrical  or  spherical  aystaa  of  coord iaatea  to  cartasian 

systan. 


(6)  BCS- 

(7)  SHAP8- 

(8)  LDB4B- 

(9)  JAC0BI- 

(10)  TI1AM- 

(11)  XXBBXL- 

(12)  mot* 

(13)  ms- 

OA)  fHAIB- 


Input,  process  aad  output  of  tha  boundary  coadltions. 

Calculates  shape  fuac tlons  aad  their  da  rivet ires. 

Calculates  linear  shapa  functions  aad  thair  derivatives. 

Evaluates  aa  Jacobian  and  tha  eoapoaaats  of  tha  unit  vector  aoraal. 

Calculates  aa  Jacobian  for  tha  case  abara  tha  first  arguasat  of  aa 
karaals  ia  a  node  of  the  alaaaat. 

Calculates  kernels  0  aad  T  . 

Carries  out  tha  auaarlcal  integration  of  kernel  products  to  fora  the 
■atria  aad  second  asufcer  uhara  tha  first  arguasat  is  not  a  coda  of  tha 
alaaaat. 

Carries  out  tha  la  t  egret  ion  of  Oe  kernel  products  uhara  tha  first  argu- 
asat  is  a  node  of  tha  alaaaat. 

Applias  tha  boundary  eoadltloas  to  feta  the  astria  aad  the  second  aaafeer. 


s 

s 
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(13)  MOT-  Mm  the  linear  nitlm  by  the  US  fewlatlM.  Ibis  miisa  truti 
eu  M  earn  only. 

(It)  OCtrot-  Outpeta  tba  Mil  UifUcmoti  tad  th«  tractlou. 


(It)  OCirot-  (tatputa  tba  Mil  dlsplaeauuts  mi  tba  tractlou. 

(17)  no-  CalenXataa  atraaau  at  aarfaea  Mai. 

(It)  CTtW-  Traufem  atraaau  la  tba  laeal  ayatu  ef  coordinates  to  tba  global  spate 

(19)  matt-  Caleolatu  tba  yriaalyal  atraaaaa  at  tba  eodas. 

(20)  ICO  I-  laata  tba  eeordlutu  of  tba  laearul  points. 


(21)  RUO-  Calealatu  tba  atrala 


af  tba  alutle  My. 


i-  Calealatu  tba  eqolaalmt  nodal  farcaa  af  tba  prescribed  ualfora  atruaaa. 
i-  Calealatu  M  utyata  tba  dlaplicaasats  and  atraaau  at  latarnal  points. 


Calealatu  tba  karaela  far  laearul  paints . 


thia  aactloa  taaarlboa  tba  operating  lu tractlou  far  tba  MX3D  prograa.  Buerlp- 
tlou  of  tba  ramlta  graiaeaf  for  a  seccaaafully  unrated  problem  ara  described  la 
SppeMix  D. 


for  all  BU  uthods  udelliag  af  yrablom  aloaya  atarta  by  disc  retiring  tba  outalda 
mrfau  ableb  la  boon  u  tba  boundary  aoatoar  of  tba  douia.  The  beoadary  contour  of  a 
tbraa  dluaolnaal  body  la  rayraautad  by  aarfaea  clients.  Tba  type  of  elauet  uad  la 
tba  pragma  at  praaaat  la  roatrictad  to  u  al^t-aoda  laoparautrlc  paadri  la  tarsi  aarfaea 
alauat.  u  abau  la  Mg  1.  Tba  parautrla  repraaaatatloa  ef  tba  alauata  enable  eorra- 
taraa  to  be  node  Had  abate  the  gaautry  ef  tba  aarfaea  la  defined  la  term  ef  guadratlc 
ah  ape  femtleu  of  tba  latrine  ie  aaacdlaatai  at  the  alauat. 


The  udellleg  teabelgu  la  uaally  staple,  bet  a  daelalon  mat  be  ude  on  boo  easy 
aarfaea  elaueta  are  uadad  to  repraaut  aay  part  lea  lar  prablm  to  tba  deeirad  accuracy. 
Thia  la  aeawlly  d  lata  ted  by  the  aa^lwity  ef  tba  atroeeaml  abaft  to  be  node  llad.  Tba 
aeeuracy  of  the  reulta  teada  to  luraoa  u  the  MM  ef  alauata  oaad  la  tba  udellleg 
ef  tba  preblaa  luroaau.  bet  tba  aaepetlag  tlu  and  out  lacrouaa  too.  ibarefora  tba 
earn  out  a  trike  a  balama  bateau  them  faetera .  Although  tba  aodalllag  of  a  prablm 
la  relatively  a  lap  la  aad  strai^ufemerd  eeu  rulu  uat  be  observed  la  order  to  preserve 
tba  eeapatibillty  ef  the  alauat  aodu  be  too  an  alauata. 

Then  an  uly  a  fee  telu  ebiab  aaed  to  be  oboonod  ehu  aodalllag  tba  aarfaea  ef 
a  fhree-dluuiaeal  body  and  then  an  daacrlbad  bo  loot 


(a)  ten  »t  l>lUt7  ef  a  leant  aadu  for  all  the  elaeeata  out  be  preserved  i 

(b)  a  aaMerlag  aaaeutlu  nut  be  observed  ebu  auMeriag  tba  surface  alaeanta. 
Tba  olfht-aodo  guadrllataral  a  leant  uat  be  defined  by  a  clockviaa  traaaraa  around 
tba  alauu  u  vieoed  u  tba  on  tilde  surface  i 


4 

Cc)  tha  eight-ooda  quadrilateral  alaaaat  mat  kna  four  aUaat  although  the  tUai 
mad  aat  to  ba  equal  la  length,  they  aaad  to  ba  eomtlauoua  aad  arneoth  ao  that  tbay 
tan  ba  fitted  accurately  with  tha  quadratic  ahapa  function. 

3.2  hlamat_aubdlvi!ioa 

day  eight-node  quadrilateral  aleaart  my  ba  tubdivided  into  a  nunbar  of  amXUr 
alamata.  far  awala,  a  largo  (quart  alaaaat  can  ba  iplit  late  two  halve!  by  drawing  a 
line  through  tha  Biddle  aa  ahem  la  fig  2,  but  la  delag  ae  we  aaad  to  latreduee  a  total 
at  five  extra  nodal  poiata  la  order  to  reatore  tha  compatibility  of  two  aight-noda 
alamata  aa  ahom  la  fig  3.  further  nodaa  my  have  to  ba  introduced  la  adjacent  elaeenta. 
for  waagle,  la  tha  aedalllag  of  a  thraa-diaaaaloaal  body,  each  aa  tha  cube  la  fig  4a, 
the  aebdlviaton  of  eoa  of  tha  element!  on  one  aide  of  tha  cube  eauaea  lacomatibillty  of 
nodal  pointt  with  tha  elaaeata  adjacent  to  it.  further  additional  nodal  point a,  in  all 
a  total  of  12  ai  ahom  la  fig  4b,  mat  ba  introduced  la  order  to  reatore  tha  compatibility 
of  nodaa.  Generally,  alaaaat  aubdivialoa  or  a  refined  aaah  will  imrovo  tha  accuracy  of 
the  reoulta  but  tha  procaea  of  elaaeat  aubdivialoa  alao  rapidly  lacreaeea  tha  amber  of 
elaaanta  for  a  three-dlaanaional  body.  Tha  extra  elaaeata  created  iacraaae  tha  aaouat 
of  input  data  aad  tha  coaputational  tine  required,  however  careful  aaah  deaign  aad  tha 
uaa  of  arbitrarily  ahapad  eight-node  elaaeata  caa  iarore  the  modelling  efficiency, 
further  laprnvumnr  la  tha  formlatioaa  of  the  ahapa  function  would  allow  triangular 
alaaaara  to  ba  aaad  and  would  provide  avea  greater  efflclaacy  la  term  of  aaah  deaign s  a 
apaclal  algorleba  la  aaadad  for  thia  aad  it  ia  being  developed  jointly  with  Southaaptoa 
del  vanity. 

3.3  Daacrlotloa  of  tha  input  data 

fur  each  problem,  the  program  raada  ia  a  earner  of  aata  af  data  ia  tha  order  ahom 
below,  tha  phyelcal  unite  for  which  mac  ba  a  eelf-cauaieteat  aat.  Apart  from  tha  TRU 
aad  tha  apecificatlaa  of  tha  coordlaata  eyatam,  all  input  data  aay  ba  catered  la  the  fora 
of  free- format.  Tha  boat  may  to  aaa  haw  thia  worka  la  to  follow  on  wangle  which  may  be 
found  la  auction  3. 1  aad  tha  corraapoadlag  input  data  file  la  Hated  la  Appendix  A  with 
deecrlptieaa.  It  aheuld  be  not ad  that  act  all  iapue  data  da flood  below  aay  ba  praaaat 
tor  a  given  problem. 

4.  ftoblam  title 

**  hater  problem  title!  it  aaa  ba  N  character!  long. 

I.  Quadrature  formic  data 

**  (near  tha  order  of  Gauaalaa  quadrature  forme  la  need  far  the  numerical  integration  over 
quadrilateral  elaaeata  aad  tha  aerreep ending  Cancel an  abaciaaaa  aad  coafficlauta. 

(hotel  Tha  Gauaalaa  abaciaaaa  aad  ita  coafficlauta  arc  uaad  la  tha  Gauaaiae  quadrature 
foraulaa  to  avaluata  varlaua  lategrala  over  alaaaata  aad  ealla  la  III.  A  abort  aumaary 
af  thane  value  a  la  Hated  ia  Appendix  t  or  a  coaprehaaaive  aat  af  the  coafficlauta 
my  ba  fewad  ia  laf  10.  Aa  a  guide  Haa  fourth-  ar  aixth-order  of  Gauaalaa  quadrature 
la  cuff iclaut  far  mat  aaaaa.  If  la  doubt,  tha  aafeat  way  to  gat  information  an  the 
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stability  of  the  results  it  to  sake  imnl  run*  with  different  order*  of  Gontiiaa 
fHintttti) 

C.  M*t«ri*l  orooortio* 

*e  Inter  the  Tounf'i  modulus  end  Poitaon'e  ratio  of  the  elastic  body. 

D'.  Ceonstrlc  nesh  date 

ee  Inter  the  total  nmdiar  of  surface  nodal  poinea,  the  total  maber  of  aurfaee  elaaenta, 
the  wash  data  output  control  parameter  (is  •  or  t)  where, 

t  t  aaah  data  output  suppressed 
1  t  aaah  data  printed  out 

and  the  aato  physical  diatnaien  of  tha  doaain  (diatance  between  the  two  furtheat 
apart  boundary  elaaenta) »  thla  value  need  only  be  eatiaated  to  e  sound  figure. 

(Rote:  Problem  requiring  total  nuaber  of  nodal  point*,  RRP<332  and  total  nuaber  of 
alaaants,  HL<t1f  can  be  treated  without  g««g  tha  diaenaiona  of  the  array*  in  the 
program.) 

I.  Rodal  coordinate*  input 

a  inter  the  type  of  system  of  coordinate*  used, 

SYS-CU  -  Cartesian  coordinate  systea 
m-CTL  -  Cylindrical  coordinate  systea 
(TS-tPI  -  Spherical  coordinate  syataa 

ee  inter  the  nodal  point  maters  with  tha  corraspoodiag  coordinates  in  accordance  with 
the  following  tablet 


*T»  - 

X(I> 

TCI) 

Cartesian 

B 

B 

Cylindrical  polar 

B 

B 

Spherical 

B 

B 

Bj 

where  the  angles  •  and  d  are  the  usual  coordinate  angles,  and  are  nea cured  in 
degreea.  In  this  Neaarandua  ananples  are  given  using  Cartesian  coordinates  only. 


f.  lleneat  data 

ee  inter  tha  element  nuaber*  cad  tha  e tenant  node  airttr*  of  the  right-node  quadrilateral 
a  tenant.  (The  coordinates  corresponding  to  the  node  nuubers  nw  previously  defined 
in  I  above.) 


6.  Reabet  of  Internal  solute 

ee  Sneer  the  total  nuaber  of  internal  point  solutions  required!  these  usually  are  the 
point*  of  interest,  inaide  the  boundary  of  the  body,  idler*  stresses,  diaplecensnta 
end  traction*  are  required.  If  no  internal  solution  is  required  then  sens  aist  be 


2 


•stand  lucud  and  skip  eo  B.  For  non-saro  yalues,  the  coordinate*  of  aaeh  interior 
point  an  to  ba  antarad.  Tha  format  of  tha  eoordiaata  17* tan  tor  tha  istarsal  point* 
■Mt  ba  tha  *aa*  a*  for  tha  aurfaca  point*. 

I.  Boundary  condition*  data 

For  data  input  purpoia,  tha  boundary  condition*  an  of  thraa  type*;  diaplaeanant 
constraint*,  neda  constraint*  and  applied  oniforn  itntaa*  oyar  certain  alaaant*. 

**  Bntar  tha  total  naahar  of  alaaant*  with  diaplaeanant  banadary  condition*. 

(Rotas  mrfnim  naahar  of  diaplacaaiat  banadary  condition*  w*t  ba  las*  than  110  if 
tha  dimensions  of  tha  array*  in  tha  pro  (ran  an  not  to  ba  chanted.) 

*»  tutor  eba  total  naahar  of  point  eanatninta  nhara  tractions  an  aaro. 

(Rotas  Point  constraint*  arc  often  nsad  in  IZX3S  to  fix  tha  rigid  body  notion  of  a 
body  la  a  chosen  direction  a*  in  tha  axanpla  in  taction  5.2.  lonayar,  tha  naar  should 
assure  ahan  cbooaiag  tha  poiats  of  constraint  for  this  purpose  that  they  will  not 
causa  any  undesirable  affects  to  tha  final  calculated  solutions.  Tha  number 

of  point  constraints  nust  ba  lass  than  30  if  tha  diaaasions  of  tha  anayt  in  tha 
prograa  are  not  to  ba  changed.) 

a*  later  tha  total  nnabar  of  sate  of  applied  oniforn  atraasas. 

If  than  an  one  or  non  alaaant*  con* trained  than  for  each  alaaant  inyolwad: 

•a  Inter  tha  constrained  alaaant  nuahar,  the  direction  ia  shich  the  alaaant  is 
constrained  (1,  2  or  3),  utsara 

1  s  proscribed  displacanmt  coup  coant*  in  tha  X -direction 

o 

2  1  prescribed  displacanmt  conponant*  ia  the  T -direction 

3  s  praacribad  diaplaeanant  conponant*  ia  tha  Z-dinctlon 

and  enter  the  praacribad  displacanant  conpooaats. 

(Rota:  The  ordering  of  tha  alaaants  is  srbltnry.) 

If  than  an  one  or  non  node  constraints,  for  each  nodal 

a*  later  tha  constrained  nodal  point  nuahar  and  enter  tha  direction  of  constraint 
(1,  2  or  3),  where 

1  1  praacribad  diaplaeanant  conponant*  is  tha  l-diraction 

2  1  praacribad  displacanant  conponant*  ia  tha  T -direction 

3  1  proscribed  displacanant  conponant*  in  the  Z-diraction 

(Rota:  The  ordering  of  tha  nodes  is  arbitrary.) 

If  thara  are  one  or  nora  sets  of  uaifota  applied  strasMi,  for  each  sat: 

••  Inter  tha  total  nu^ar  of  alaant*  sub j  acted  to  tha  uniform  applied  stresses 

and  specify  tha  nognituda  of  tha  stress  field*  in  Ps,  Fyy,  Fss,  Fry,  Put  and  Fyt. 
(Rota  1  A  combination  of  these  sera**  field*  nay  ba  specified  so  that  tha  slnulaead 
load  can  ba  characterised . ) 
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**  bear  Ckt  total  number  of  alaaanti  which  ara  got  orthogonal  to,  but  ahara  nodes 
with,  thoaa  alaaanti  aetad  on  by  thaaa  atraaaaa. 

(Motel  Tba  and—  nuabar  of  alananta  subjected  to  tha  uniform  applied  atraaaaa 
scat  ba  laaa  than  M#  tf  tha  dimensions  of  tha  arrays  ta  tha  program  ara  not  to  ba 
ehangad.) 

**  Intar  tha  alanant  number  of  tha  alananta  that  ara  subjected  to  tha  above  aat  of 
uniform  loads. 

up  of  m« 

3.*  Mssh-noaaration  program  and  graphic  plotting  co^utar  orogran 

llthoapi  tha  input  ia  auch  ai^lat  for  tha  boundary  alanant  nathod  than  for  tha 
finita  alanant  nathod,  neat  of  thia  advantage  would  ba  loat  if  tha  coordlaataa  of  aach 
naah-point  and  tha  connectivity  af  aach  alanant  bad  to  ba  input  aaparataly.  Tharafora, 
in  oxdar  to  naka  tha  precaaa  of  c ranting  a  data  fila  aa  (lap la  aa  poaaibla,  a  threa- 
diaanaional  waih-gana ration  program  haa  baan  developed  and  iaatallad  ia  tha  Materials  and 
Structural  dapartnantal  UX  eonputar.  The  program  not  only  lavas  tina  but  alao  alinlnatas 
nasty  of  tha  human  arrora  which  nay  occur  during  manual  preparation  of  tha  coordlaataa  and 
alanant  topologlaa.  lha  nesh-gaaeratioa  program  ia  callad  MESH3D.  It  la  based  on  auto- 
natic  alanant  aubdlwialon  of  af  aw  largo  blocks  which  ara  daflnad  as  Input  datai  aach 
block  is  than  subdivided  into  alananta  according  to  tha  subdivisions  which  ara  provided 
as  input  data.  1CSH3D  will  calculate  all  the  nodal  coordlaataa  for  all  tha  subdivided 
alananta.  Tha  nodal  point  nuabar  log  convention  (clockwise  or  anticlockwise),  generated 
for  tha  mash,  will  ba  tha  same  aa  the  input  data  nunbsring  convention.  At  present,  tha 
naahes  calculated  by  tha  MZSH3D  are  ia  Cartesian  coordlaataa  only. 

After  a  nash  has  baan  gaaaratad  by  tha  MBSI3D,  it  can  ba  visually  checked  by  plot¬ 
ting  out  tha  mash  on  a  graphics  display  screen  or  on  paper.  The  plotting  of  tha  nash  is 
aa  essential  part  of  the  us ah-gaae ration  procedure  because  the  only  way  to  check  tha 
topology  of  tha  nash  quickly  and  accurately  is  by  visual  inspection  of  a  three-dimensional 
view  of  tha  nesht  tha  large  volume  and  tha  thrae-d loans ional  complexity  of  tha  nash  data 
often  nake  it  virtually  impossible  to  check  numerically. 

A  full  report*  and  a  coaprahaasiva  user’s  guide  for  both  tha  nesh-geaeratloe  program 
and  tha  graphic  plotting  program  have  been  prepared  at  tha  earns  time  as  the  present  Techni¬ 
cal  Memorandum.  Details  of  tha  basis  of  tha  mesb-geaaratieo  procedure,  subdivision  pro¬ 
cess  and  tha  connectivity  .of  individual  blocks  are  explained  with  anemias.  In  that 
report,  tha  operations  of  tha  graphic  plotting  program  sad  data  handling  are  also  explained 
with  an tuples. 

A  ggglCAL  ITALPAXIOW  OP  Tto  HUM  mBWTT  FACIOt 

The  behaviour  of  a  fatigue  crack  ia  a  given  material  is  controlled  by  tba  parameter 
callad  tha  strata  intensity  factor,  I ,  which  is  a  measure  of  tha  severity  of  the  stress 
field  at  tha  crack  tip.  Tha  strass  intensity  factor  Ij  under  Mode  X  conditions  is  a 
function  of  tha  loading  and  the  geoaetry  af  tha  crack  and  tha  body.  Tha  numerical  evalu¬ 
ation  of  crack  front  stress  intensity  factors  using  tha  III  nathod  involves  two  phases t 
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(a)  carrying  out  a  stress  analysis  of  a  cracked  body  using  a  boundary  ala 


ash. 


(b)  development  of  techniques  for  the  determination  of  tha  serosa  intensity  factors 
after  all  eha  displacements  and  erections  era  calculated  for  all  tha  alasmnts. 


It  should  bo  noted  that  la  eha  prosenc  program  (SZE3D)  only  those  crack  problems 
symmetric  about  tha  plana  of  tha  crack  and  under  Mode  I  loading  can  be  analysed;  advant¬ 
age  is  taken  of  tha  syaaatry  about  ehe  crack  plana,  thus  reducing  tha  physical  problem 
siaa  to  ba  analysed.  In  eha  modelling  of  the  physical  problem,  the  crack  plana  is 
represented  as  a  boundary  with  appropriate  boundary  eomditioaa.  The  crack  face  is 
normally  free  of  constrainta,  but  stresa  constraints  say  ba  applied  if  required.  Zero 
normal  displaeomsnts  are  prescribed  on  tha  plane  ahead  of  the  crack. 

Various  techniques  for  determining  Kj  from  the  boundary  integral  results  can  be 
used  but  is  the  present  study  only  the. ’displacement  method’  will  be  need;  it  is  explained 
in  the  following  section. 


4.1  The  displacement  method 

for  a  general  three-dimensional  crack,  the  stress  intensity  factor  will  vary  with 
the  position  along  tha  crack  front,  for  an  infinite  three-dimensional  body  tha  relation¬ 
ship  between  I  and  displacement  1m  a  plana  perpendicular  to  the  crack  front  is  the  same 
as  in  the  two-dimensional  plana  strain  case.  It  ia  therefore  assumed  that  the  relation¬ 
ship  between  K  and  displacement  in  the  plana  perpendicular  to  the  crack  front  is  given 
by  tha  two-dimensional  plana  strain  equations  if  the  plane  is  within  the  body.  For  a 
surface  plane  ie  ia  assumed  that  the  relationship  is  given  by  tho  plana  stress  equations. 
Thaaa  procedures  have  been  discussed  by  Cruse*. 


The  two-dimensional  displacement  field 
r  from  the  crack  tip  under  Mode  I  condition 


where  suffix  1  and  2  refer  to  tho  direction 
modulus,  and 


[U|,u2]  ,  as  shown  in  Fig  5,  at  a  distance 
it  defined  aa*i 


-  1  ♦  2  sin 


l] 


(1) 


♦1-2  COS*  yjj 

x  and  y  respectively,  u  is  the  shear 


(3  -  *v> 


for  plane  serein 
for  plana  stress 


(2) 


where  v  is  Foisson's  ratio.  For  points  along  a  radial  line  normal  the  crack  front  and 
in  tha  plane  of  the  crack,  the  displacsnene  normal  to  the  crack  face  (is  1  •  «  in 
Fig  ))  1st 
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where  r  is  tha  distance  from  th*  crack  front,  and  t'  is  eh*  aacarial  stiffness,  equal 
to  th*  Young's  aodulua,  K  ,  tor  plan*  ier«aa  condition*  and  equal  to  t/(t-v  )  for 
plan*  aerain. 

Zf  Eha  computed  value  of  u_  ae  radial  diieane*  r  fron  th*  crack  front  it  v  , 

*  * 

th«  corresponding  mlus  of  9  toy  Kj  ,  is  obtninsd  from  (3),  w 

^  ■ 

Th*  above  aquation  (3)  ia  strictly  valid  only  in  tha  liait  when  r  approaches 
**ro.  Thus,  values  of  Kg  ar*  obtained  using  (4)  fron  th*  computed  displaceaents  V 
at  a  series  of  nodal  points  along  a  radial  line  on  th*  crack  face  and  extrapolated  to 
r  •  0  to  yield  tha  crack  front  stress  intensity  factor  Kg  •  The  cosputed  values  T 
in  regions  very  dose  to  eha  crack  front  are,  however,  less  reliable  and  it  is  therefore 
necessary  to  extrapolate  fron  relatively  large  values  of  r  •  This  inaccuracy  is  clearly 
shown  in  the  worked  exasple  described  in  section  5.2  (Pig  13)  where  K^  is  plotted 
against  (r/C)  ,  and  where  C  ia  the  crack  length.  Th*  constant  slope  portion  of  th* 
curs*  is  extrapolated  to  r  ■  0  . 

The  extrapolation  netted  or  the  displacanent  netted  above  does  not  contain  any 
special  provision  in  the  nunsrical  fomilation  for  the  r*  dependence  of  th*  displace- 
nent  or  the  r  ^  singularity  in  traction  at  th*  crack  front,  although  the  boundary 
integral  aquation  nathod  gives  good  resolution  of  th*  solution  paraaeters,  a  relatively 
fine  grid  around  tha  crack  front  ia  still  necessary  in  order  to  obtain  reasonably  accurate 
values  of  Kg  by  these  extrapolation  netted*.  To  overcons  this  nejor  disadvantage,  it  is 
possible  for  the  users  to  enploy  th*  'quarter  point  aleanmts'  at  th*  crack  front  in  th* 
nsah  design.  The  basic  concept  underlying  the  quarter  point  elenent  is  that  th*  require- 
nent  for  th*  r^  dependence  of  th*  displacanent  near  tha  crack  front  in  the  8-node 
quadrilateral  elenents  can  be  achieved  by  placing  th*  nid-sid*  nodes  of  th*  elenent* 
either  side  of  tha  crack  front  at  th*  quarter  point  instead  of  the  mid-point.  Th*  pos¬ 
itioning  of  nodes  for  quarter  point  elenents  is  Illustrated  diagrannatically  in  Fig  4. 

It  is  inportant  to  not*  that  this  procedure  does  not  require  special  shape  functions  in 
the  iaoparaaetrie  coordinate  specs  and  therefor*  does  not  require  sny  different  conpueer 
subprogram.  The  resolution  and  accuracy  of  eh*  solution  parameters  nesr  eh*  crack  tip  can 
be  improved  further  still  by  eh*  application  of  th*  'Singularity  llenants'  ahead  of  th* 
crack  front  ia  the  nunsrical  formulae ion  of  eh*  fraceur*  problem. 

4.2  The  'Ilnsularlty  llenants1 

Th*  special  'Singularity  Ilsments',  which  exhibit  th*  stress  (strain)  and  traction 
singularities  occurring  at  the  craek  front,  have  been  used  la  tvo-diaentioual  problems. 

They  require  less  mssh  refinement  around  th*  crack-tip  region  than  th*  extrapolation  { 

method,  described' earlier,  to  obtain  accurate  values  of  th*  stress  intensity  factors. 

? 

i 
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This  concept  hot  been  introduced  fro*  the  finite  element  method  end  extended  to  the 
bounder?  integral  equation  Method  with  the  it opereaetric  forauUtlon.  The  beeic  concept 
underlying  the  eintulerity  it  that  there  ie  a  requirewnt  for  an  r“*  eingulerity  in 
traction*  which  can  be  achieved  by  imposing  a  shape  function  in  the  numerical  formulation 
to  characterise  thi*  r"*  behaviour  near  the  creek  tip.  The  special  formulation  needa 
to  be  contained  in  a  computer  subprogram.  However,  singularity  elements  in  tractions 
have  been  fowlated  in  two  dimensions  only^  and  no  three-dimensional  procedures  are  yet 
available. 

5  TEST  CASES 

This  section  contains  two  worked  examples  which  will  be  used  to  illustrate  the 
three-dimensional  boundary  element  program  (BIE3D)  and  how  it  can  be  used  to  calculate 
stresa  oonoantration  factors  and  stress  intensity  factors. 

In  order  to  show  the  method  of  data  preparation  and  the  format  layout  of  an  input 
data  file,  the  input  data  file  of  exa^le  (1)  ia  shown  in  Appendix  A  and  the  correspond¬ 
ing  output  data  file  is  shown  in  Appendix  B. 

The  BIE  solutions  are  compared  with  known  accurate  solutions. 

5.1  Determination  of  stress  oonatntratim  factor* 

A  thick  flat  plate  with  a  circular  hole  normal  to  the  plate  surface,  as  shown  in 
Tig  7,  is  a  classical  stress  omomtratien  problem  and  well  documented  in  the  engineering 
field.  Many  accurate  numerical  solutions  for  this  specimen  are  readily  available  for 
comparison  of  results.  Therefore  it  will  be  used  here  as  a  test  case  for  the  BIE3D 
computer  program. 

Tor  simplicity,  only  a  quarter  of  the  flat  plate  will  be  modelled  by  making  full 
usa  of  the  sysetry  about  the  x-  end  y-axes .  The  three-dimensional  mesh  representation 
of  the  quarter  plate  consists  of  a  total  of  52  quadrilateral  surface  elements  and  a  total 
of  158  nodal  points.  The  numbering  system  of  the  mesh  for  the  surface  elements  and  the 
nodal  points  is  shown  in  Tigs  8  and  9  respectively.  The  corresponding  input  data  file 
(TP. DAI)  is  listed  in  Appendix  A.  This  input  data  file  may  be  either  typed  in  by  the 
user  following  the  steps  described  in  section  3.3  or  it  can  be  generated  using  the  mash- 
generation  program  ME5H3D  a*  briefly  mentioned  ia  section  3. A;  as  a  demonstration,  the 
corresponding  MESH3D  input  data  file  is  shown  her*  ia  Appendix  C.  By  eonparing  the  two 
files  in  Appendix  A  and  Appendix  C  it  can  be  seen  that  a  large  saving  ia  terms  of  data 
preparation  is  possible.  A  user's  guide  for  the  preparation  of  a  MESH3D  input  data  file 
has  been  prepared  la  the  form  of  a  Technical  Memorandum3. 

finally,  the  coaamnds  required  to  submit  the  file  (TP. DAT)  for  the  BIE3D  execution 
run  are  at  shown: 

tASSXGH  TP. DAT  F0EM5  (assigned  to  logical  unit  5) 

$ ASSIGN  TP. EES  FOB986  •  (assigned  to  logical  unit  6) 

$ BUN  BXE3D 

where  TP. EES  ie  the  output  file  of  the  results,  as  shown  in  Appendix  B;  it  contains  the 
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displacements,  the  tractions  and  the  atraaaaa  for  all  cha  nodes  and  all  the  intanal 
points  (if  specified)  of  tha  quart  a  c  flat  plata.  In  thia  problaa,  tha  — atraaaaa 
ara  tanalla  in  tha  direction  of  tha  applied  atraaa  and  occur  at  tha  aidaa  of  tha  hole  on 
tha  net  aection.  The  atraaa  oonemtration  factor  can  be  defined  aa  k(  ■  o^/o  . 

In  order  to  demonstrate  tha  variation  of  atraaaaa  through  tha  thickness  of  tha 
tpaciaas  the  atraaaaa  along  one  aide  of  tha  hole  ara  plotted  froa  ana  aide  of  tha  surface 
to  the  other  as  shown  in  Pig  10.  the  maxima  atraaa  oonoentration  factor  occurs  at  tha 
aid-plane  having  a  value  of  3.26,  decreasing ‘to  3.01  at  tha  free  surface.  Also  in  Pig  10 
ara  the  solutions  obtained  bp  Sternberg  and  Sadowaky*  using  aa  analytical  approach  for  an 
infinite  plata  having  tha  same  ratio,  T/D,  of  plata  thickness  to  hole  diameter.  The 
corresponding  analytical  solutions  art  3.11  and  2.78  respectively.  Because  of  tha  effect 
of  Cha  finite  width  of  Cha  plate,  tha  8X8  solution  should  be  reduced  by  approximately  IX 
for  coaqpariaon  with  an  Infinite  plata.  Bowevar,  tha  BIE  solutions  remain  approximately 
3-4X  higher  than  those  obtained  analytically!  this  is  satisfactory  agreement  in  view  of 
the  relatively  coarse  boundary  ms  ah  used. 

3.2  Determination  of  atraaa  intensity  factors 

The  compact  tenaion  specimen  used  for  fracture  toughness  tasting  (ASTM  8399-74) 
was  represented  by  a  Single-Edge-Crack  (SEC)  specimen,  aa  shown  in  Pig  11,  with  crack 
length,  C,  width,  W,  thickness,  T,  and  length,  2B,  such  that  CH  -  T/V  -  0.3  and  H/W  •  0.6 
Poisson's  ratio  was  taken  to  be  0.3.  Taking  advantage  of  the  two-plane  symaetry,  only 
one-quarter  of  the  specimen  was  modelled.  The  boundary  mesh  idealisation  is  shown  in 
Pig  12.  It  consists  of  64  surface  elements  and  194  nodes.  The  numbering  systems  of  the 
meshes  for  tha  surface  elesstnta  and  the  nodal  points  are  shown  in  Pigs  12  and  13  respec¬ 
tively.  In  order  to  aisailate  the  coop  act  tension  specimen  load,  a  uniform  shear  stress, 
o^  ,  was  applied  to  tha  x  -  8  plane,  as  shown  in  Pig  It,  to  produce  a  total  vertical 
load  in  the  y-directioni  where  the  uniform  shear  stress  used  in  thia  present  ezaaple  was 
unity. 

The  stress  intensity  factors  along  the  crack  front  were  obtained  by  the  extrapol¬ 
ation  method  with  tha  use  of  crack  face  displacements  as  described  in  section  4.1. 

Por  the  present  exaaple,  tha  crack  face  displacements  at  tha  tree  surface  plane,  mid- 
plane  and  symmetry  plane  ware  used  in  order  to  calculate  the  stress  intensity  factors 
along  the  crack  front.  The  results  are  presented  in  graphical  form  in  Pig  14  wbera  the 
normalised  stress  intensity  factors  (E^/fy(C)^/TWJ)  are  plotted  against  (r/C)  and  the 
linear  portion  of  the  eurve  ia  then  extrapolated  to  the  erack  eip  and  yields  tha  value  of 
Kj  at  r  •  •  .  It  should  be  noted  that  the  reaults  shewn  in  Pig  14  were  obtained  with 
the  conventional  iaoparanetrie  elements.  Plans  strain  conditions  wars  assumed  to  prevail 
at  all  points  along  the  crack  front,  except  where  it  intersects  ehe  free  surface  and 
there  plane  stress  assertion*  were  nada. 

The  variation  of  the  normalised  stress  intensity  factors  along  tha  eraek  front  are 
plotted  in  Pig  tS.  Also  shown  ta  Pig  IS  is  tha  two-dimensional  plane  strain  solution 
tor  the  actual  teat  specimen  geometry  obtained  by  Newman  (1974)*  using  tha  boundary 
collocation  method,  which  is  13.68.  Tha  HE  solution  obtained  ia  this  case  is  13.36, 
which  is  about  2Z  lower  at  the  aywatry  plane  than  Newman’s  result. 
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ImUit  ftoklM  «w  Htlyiid  tulng  the  ua  specimen.  la  this  east  a  imifon 
tensile  s trass,  ,  ass  applied  at  the  7*1  plans,  as  shown  in  Fig  It.  Th«  strsss 
intensity  factors  wars  again  calculated  using  tha  dlsplacaaant  aathod  and  thair  variation 
across  tha  thichaass  of  tha  spsffaan  la  shown  la  Fig  16.  Ihs  stress  intensity  factors 
wars  aatwallsad  by  OqIvC]^  la  this  case,  and  at  tha  syaastry  plana  tha  mine  obtained 
was  2.65,  while  that  at  tha  free  surface  was  2.60.  Again,  pleas  strain  conditions  wars 
assnaad  to  preoall  at  all  nodal  points  along  tha  crack  length,  eaeept  at  tha  free  surface. 

Bo  results  far  tha  aaaa  gaoostry  aad  loading  are  available  la  tha  literature  for 

comparison,  however,  a  two-dimensional  plana  strain  solution  of  the  noraalisad  K.  of 

1  £ 

2.60  have  bean  calculated  for  a  speciaea  with  tha  saaa  crack  length  by  Bowie  and  leal 
using  tha  coaforaal  napping  technique. 

6  PUCT8SI0M  MB  CBTOBMOW 

la  tha  present  Mameraadum  a  casp liter  progran  based  on  tha  Boundary  Integral 
Equation  aathod  has  boon  described.  The  essential  parts  af  tha  operating  procedure  of 
the  progran  have  bean  demonstrated  with  enables.  It  has  been  shown  that  with  care  the 
conputer  progran  can  be  used  as  a  numerical  tool  to  calculate  tha  stress  intensity 
factors  for  a  cracked  body,  the  method  illustrated  in  tha  worked  anaaples  for  the  cal¬ 
culation  of  the  stress  intensity  factors  la  based  on  tha  displacement  extrapolation 
aathod.  Mathematically,  the  extrapolation  aathod  Is  based  on  the  two-dimensional  stress 
fields  In  a  two-dimensional  plana  body,  therefore  it  is  net  an  exact  analogy  in  the 
thxee-dimensloual  cases.  However  for  an  Infinite  three-dimensional  body  the  relation¬ 
ship  between  K  and  displacement  in' a  plana  perpendicular  to  tha  crack  front  is  the  earns 
as  in  the  two-dimensional  plana  strain  case.  On  the  other  hand,  tha  surface  of  a  three- 
dimensional  body  cannot  be  represented  by  the  two-dimensional  plane  stress  equations. 
Therefore  bettor  methods  are  needed  so  that  more  reliable  stress  intensity  factor 
solutions  can  be  obtained  in  tha  future.  Currently,  there  are  two  methods  being 
developed  jointly  in  IAE  Farnborough  aad  the  Mechanical  Engineering  Department, 
Southaapton  University: 

(t)  Subtraction  aathod1 1 ,  and 

(2)  Buackner  singularity  method'*. 

A  three-dimensional  mssh-gensraeion  program3  has  been  developed  as  a  tool  to  pre¬ 
pare  Input  data  files  for  the  BXUS  computer  program.  When  the  mash-generation  progran 
is  coupled  with  the  BXZ3D  program,  the  user  can  solve  problems  with  greater  efficiency. 

In  order  to  obtain  a  visual  check  of  tha  three-dimensional  view  of  tha  mesh,  a  three- 
dimensional  graphic  plotting  program  has  also  been  developed.  Theqe  computer  programs 
provide  a  useful  part  of  the  software  support  for  the  future  development  of  any  numerical 
stress  analysis  programs  to  obtain  stress  intensity  factors  for  cracks  in  three- 
dimensional  bodies. 
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-0.96961*00 

hbal  buflacbbhti 

■OH  OF* 

1  .00001*00  ..12001-02  .91906-04 

9  .00001*00  -.111*1-02  .1 90*1-04 
9  .00001*00  -.  1 16*1-02  ..1964B-0* 


199  -2121B-09  99*21-03  .89T71-09 
199  .21211-09  -.99*21-03  -.69*61-09 
197  .99921-04  ..109*1-02  .64971-09 


itumh  n  m  global  cabtbbiab  stum  or  columns 
■or  non  non  nos  sxozr  non  axon 
1  0.26991*00  0.10621*01  -0.20961*00  -0.39991-01  0.21131-01  0.6*606-02 
2  0.92991*00  0.10*91*01  -0.12*71*00  0.00001*00  0.7*306-02  0.79201-02 


■OH  0  T  V 

2  .00006*00  -.11001-02  .29631-0* 

4  .00006*00  -.11671-02  .00691-06 

(  .00006*00  -.11001-02  -.29611-0* 


194  .20991*09  -99971-09  .44*91-06 
19*  .21066-09  -.99421-09  -.17211-0* 
196  .99921-04  -.10961-02  -.64221-09 


to  m  m 

0.  0.  0. 
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157  0.1252*00  0.1000*01 

158  0.12528*00  0.10008*01 


0.12458*00  0.00008*00  0.0000*00  -0.88058-02 
0.12458*00  0.00008*00  0.00008*00  0.88118-02 


macmi.  mum  as  rat  boms 


BOM 

8101 

8X02 

8X05 

1 

0.10858*01 

-0.21048*00 

0.28218*00 

2 

0.10458*01 

-0.12698*00 

0.52568*00 

5 

• 

• 

0.10518*01 

-0.15578-01 

• 

• 

0.50558*00 

ni 

0.10558*01 

• 

-0.58658-01 

-0.19898-01 

157 

0.10008*01 

0.11598*00 

0.15568*00 

158 

0.10008*01 

0.11598*00 

0.1556B*00 

RBA1B  888801  or  188  BOR 

BOS  888MB  18  ABAX.TSXS. 

3 

9 


»/ot 


#5ia!aKSS!i8»S3«:c;S«i32o*“0,'4O",,*‘',Mtf 


tw ndi»  C 

«SH3D  noga  UffPT  DMA  HU 
(•«•  •action  5.1) 


nil ii  rut  run  vita  a  cams  sou 

6 

.25861919  .46791395 

.66120940  -56076158 

.95246949  .17132449 

-.23861919  .46791393 

-.66120940  .36076158 

—93246949  .17132449 

1000.0  0.3 

32  10  8  3 

1  12345«78 

2  5  4  3  9  10  11  12  13 

3  12  11  10  14  15  16  17  18 

4  17  16  15  19  20  21  22  23 

5  7  6  5  13  12  24  25  26 

6  25  24  12  18  17  23  22  27 

7  7262528293018 

8  25  27  22  21  20  31  29  28 

9  20  19  15  14  10  32  29  31 

10  1  30  29  32  10  9  3  2 

1  10  0 

1  0  .1 
1  0  .2 
1  .4  .2 
1  .8  .2 
1  .8  .1 
1  .8  0 
1  .4  0 
.5  0  .2 
0  0  .2 

.4293  .4293  .2 
.8586  .8586  .2 
.9235  .8152  >2 
0  .5  .2 
0  1  .2 
.4  1  .2 
.8  1  .2 

.8152  .9235  .2 
0  1  .1 
0  1  0 
.4  1  0 
.8  1  0 
.8  1  .1 

.8586  .8586  .1 
.8586  .8586  0 
.9235  .8152  0 
.8152  .9235  0 
.4293  .4293  0 
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29  0  0  0 

30  .5  0  0 

31  0  .3  0 

32  0  0  .1 

1  3  4 

3.5  3  3-5 

2  3.5  3  1.5 

2  4  1 

1.5  3  3-5  2 
1 

3  4  1 

1.5  3  3-5  2 
1 

4  4  3 

1.5  3  3.5  2 

3.5  3  3.5 

5  3  1 
3*5  3  3-5 
1 

6  3  1 

3.5  3  3-5 
1 

7  1  4 
1 

1.5  3  3.5  2 

8  1  4 
1 

1.5  3  3.5  2 

9  3  1 

3.5  3  3-5 
1 

10  1  3 
1 

3.5  3  3-5 

S 
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Appendix  p 

DESCRIPTIONS  OF  THE  RESULTS  OF  A  PROBLEM 

lull  descriptions  of  as  output  file  of  a  successfully  executed  problem  are  listed 
below  la  this  aoetioa.  Each  individual  output  deecription  in  as  output  file  will  be 
explained  eo  that  the  uaera  can  follow  tha  output  daacriptiona  and  cheek  against  their 
own  filoa  if  required. 

Descriptions  produced  by  the  cooputer  prog ran  BI13D  will  be  quoted  here  aa 
< output  deecripeios>. 

Ike  raaulta  produced  for  a  eueeeaafully  executed  problem  are  printed  in  the 
following  order.  The  units  axe  those  of  the  input  data  which  aust  be  a  consistent  set. 

Problem  title  is  printed  according  to  the  particular  problem  title  supplied  as  input 
data. 

Gaussian  quadrature  formula 
A  heading: 

< specified  gauss  abscissae  abd  coefficients  for  elemental  integration  in  the  local 

COORDINATE  SYSTEM) 

is  printed  followed  by  a  headed  table  of  the  Gaussian  abscissae  and  coefficients  supplied 
as  input  data. 

Material  property  data  . 

The  Young's  modulus  and  Poisson's  ratio  of  the  elastic  body  are  printed  following 
the  aessages  <  YOUNG'S  MODULUS  ">  and  <  POISSON'S  RATIO «  >  respectively. 

Mesh  data 

If  the  aash  data  output  control  paraaster  (MOOT)  has  the  values  taro,  no  aesh  data 
are  printed.  The  value  of  this  paraaater  should  be  defined,  normally  by  being  read  in 
as  data,  in  subprogram  MESH.  Any  non-taro  value  causes  a  full  set  of  aash  data  to  be 
printed  in  the  following  order: 

(i)  A  heading  CGEOMETRIC  DATA  FOR  THE  MESH) 

(ii)  A  asssage  <  NUMBER  OF  ELEMENTS  •  >  followed  by  the  relevant  number. 

(ill)  A  aeasage  <  NUMBER  OF  NODAL  POINTS  ■  >  followed  by  tha  relevant  nuaber. 

(iv)  A  headed  table  of  node  nuaber  and  global  coordinates  for  each  node  in  turn  in 
numerical  order. 

(v)  A  heading  CELEMENT  DATA>  followed  by  a  headed  table  e.'  element:  nuaber  and  the 
nuabert  of  its  eight  nodes. 

Boundsry  conditions 

(i)  A  heading  DISPLACEMENT  BOUNDARY  CONDITIONS)  is  printed.  On  the  next  line  a  heading 
< CONSTRAINED  ELEMENTS)  is  printed,  followed  by  a  headed  table  showing  eleaene  nuaber, 
restraint  condition  type  nuaber  (I  or  2  or  3)  where 
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in  the  X-diraction 
in  the  Y-direction 

in  the  Z -direction  and  the  coiponents 

(ii)  A  heading  <TRACTION  BOUNDARY  CONDITIONS)  ia  fallowed  on  the  nest  line  by  one  or 

■ore  headed  table (a)  (hewing  the  total  nuabar  of  eleaanta  under  each  uniform  load¬ 
ing,  tha  three  atraas  c opponent a  and  the  nunbera  of  the  elcaenta  over  which  the 
atreaaaa  act.  Thi*  ia  followed  by  a  headed  table  ahoviag  the  elaaent  nuabar,  the 
node  nuaber,  condition  type  nuabar  (f  or  2)  where 

t  s  preacribed  traction  coapoaanta  in  the  X-direction 

2  i  preacribed  traction  coaponanta  in  the  T-direction 

3  :  preacribed  erection  coaponanta  ia  the  Z-direction  and  the  coaponente  of 

the  preacribed  traction*. 

Model  point  tract iona  and  diaplactaanta 

(i)  A  heading  <CALCOLAIED  VALUES  0 F  TIACTIORS>  la  printed  followed  by  a  beaded  table 
ehowing  the  node  nuaber,  condition  type  nuabar  (I  or  2)  where 

t  :  preacribed  traction  coaponanta  ia  the  X-diraction 

2  :  preacribed  traction  covenant  a  in  the  T-direction 

3  :  preacribed  traceion  coaponanta  in  the  Z-direction  and  the  coaponanta  of 

the  calculated  tractiona. 

(ii)  A  heading  <HODAL  DISPLACEMENTS)  ia  printed  followed  by  a  headed  table  ahowing  the 
node  nuaber  and  the  computed  diaplaceaant  coaponanta  for  each  node  in  turn  in 
numerical  order. 

Shifted  node* 

-t/2 

In  order  to  iapoaa  the  r  cingularity  at  the  crack  tip,  tha  aid-aide  node*  of 

tha  aleaentt  adjacent  to  the  crack  tip  are  ahifted  to  quarter  point*  neareat  to  the  tip; 
bonce  tha  Jacobian  of  the  tronaforaaeion  froa  global  to  local  coordinate  ayatea  becoae* 
aero.  Thia  caaea  tha  following  aaaaaga  to  be  printed. 

<JAC0B1AN  »  0  .  AX  MODE  •  HN  ELEMENT  HO  ■  WO 

with  MM  and  MM  replaced  by  the  relevant  nunbera.  Tha  pregraa  ia  not  halted.  However, 
the  atreaaaa  at  that  node  are  not  calculated  aince  the  (treaaaa  would  be  infinite  there. 

Modal  point  atreaaaa 

A  heading  < STRESSES  IM  THE  CLOBAL  CARTESIAN  SYSTEM  OF  COORDINATES)  ia  printed 
followed  by  a  headed  table  ahowing  the  node  nuabar  and  the  SIX  atreaa  coaponanta  for  each 
node  in  nuaarieal  order. 

Interior  point* 

If  any  internal  aolution*  are  required  a  heading  < INTERNAL  POINT  COORDINATES)  i* 
printed  followed  by  a  headed  table  ahowing  internal  point  nuaber  and  the  coordinate*  of 
the  point. 


1  s  preacribed  diaplacaaent  component 

2  >  preacribed  diaplacaaent  component 

3  s  preacribed  diaplacaaent  covponent 

of  the  preacribed  diaplacaaent. 
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Plsglacaaent  And  *tr**»«*  »t  Interior  point! 

(i)  Die  first  »ee  of  the  final  eet  of  remit*  to  be  written  are  the  etreeeee  at  internal 
point*.  A  heading  < STRESS  AX  OTESRAL  F0XHTS>  i*  followed  by  a  headed  table  showing 
the  interior  point  nuaber  and  the  three  stress  components  for  eaeh  point. 

(ii)  Finally,  a  heading  DISPLACEMENTS  AX  INTEBNAL  P0INTS>  is  followed  by  a  headed  table 
showing  the  interior  point  nuaber  and  the  displacement  components  for  each  point  in 
turn. 


Aggandj*  E 


anmUT  or  THE  C4PSSIAH 


peatoee  cotmciprrs 


Gauaaian  abaci****  t  Oauaalaa  coafficianta 


i  0.3599610435  0.6521451548 

4  i  0.6611363115  0.3478548451 

i  -0.3399810435  0.6521451548 

i  -0.6611363115  0.3478548451 


t  0.2366191860  0.4679139545 

l  0.6612093664  0.3607615730 

6  l  0.9324695142  0.1713244923 

t  -0.2386191860  0.4679139345 

i  -0.6612093864  0.3607615730 

i  -0.9324695142  0.1713244923 


i  0.1834346424 
i  0.5255324099 
«  0.7966664774 

8  <  0.9602896564 

»  -0.1834346424 
i  -0.5255324099 
i  -0.7966664774 
l  -0.9602898564 


0.3626837833 

0.3137066458 

0.2223810344 

0.1012285362 

0.3626837833 

0.3137066458 

0.2223810344 

0.1012285362 


ft/ftl  W  HI 
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Fig  8  Mesh  design  snd  elewent  numbering  system  of  the  fist  piste 
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Fig  9 


Fig.  9  Mesh  design  and  nodal  numbering  system  of  the  flat  plate 
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SHEAR  STRESS  CASE 


SSS 
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UNIFORM  TENSILE  STRESS  CASE 


CRACK 

SURFACE 


'  W 


Fig  11  Single  edge  crack  (SEC)  specimen 
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Fig  12 
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Fig  12  Mesh  design  and  element  wintering  system  of  the  compact  tension  specimen 
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Fig* 13  Mesh  design  and  nodal  numbering  system  of  the  compact  tension 
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SURFACE  PLANE 

INTCP-  11.75 


MID-PLANE 

INTCP-  13.22 


SYMMETRY  PLANE 

INTCP-  13.38 


Fig  14  Normalised  stress  Intensity  factors  (K*/[F/C/TN])  against  (r/C) 
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Fig  15  Variation  of  normalised  stress  Intensity  factor  (kT/CF/E/TV]) 
with  thickness  for  SEC  specimen 
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Fig  16  Variation  of  normalised  stress  Intensity  factor  (K*/[o/rC]) 
with  thickness  for  SEC  specimen  (uniform  tensile  load) 
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